Over recent years the advancement made in computing power has revolutionised the process of designing aerodynamic structures. Traditional trial and error methods have been replaced by automated optimisation routines. The general idea is that by parameterising a design by a series of design variables and then automating design generation and assessment, a search process can be used to find optimum designs. During such an optimisation process hundreds/thousands of designs can be produced that have to be evaluated. It is therefore important to keep the computing resources needed per design to a minimum. However, it is now also commonplace within the aerodynamics industry to employ high fidelity computational fluid dynamics (CFD), a computationally expensive method, as a design quality assessment tool.
A typical optimisation process using CFD analysis is illustrated in FIG. 1. Here, a CFD mesh representation 2 of an aerodynamic structure is generated from a corresponding Computer Aided Design (CAD) model 1. A CFD analysis is performed on the mesh at step 3 and the results are analysed at step 4 to determine whether the design has been optimised. If so, the optimal design is output from the process at step 5. If not, the design is altered at step 6. This will be explained in more detail below. A CFD analysis is then performed on the altered design at step 3 and the results are analysed again at step 4. This process is repeated until an optimal design has been found.
For aerodynamic shape optimisation a flexible parameterisation method is required to generate new aerodynamic designs at step 6. Within the field of shape optimisation the method of shape parameterisation is critical to the overall success of the optimisation process. The parameterisation technique defines the design space and hence determines the range of possible designs the optimisation process can produce. The suitability of a particular parameterisation technique is based on its ability to balance the need for a flexible representation with a manageable number of design variables. If the parameterisation is too inflexible the number and variation of the possible shapes is restricted, thereby restricting the outcome of the optimisation process. Conversely, a very flexible parameterisation is likely to require a large number of design variables. If the number of design variables is too large the design space becomes inefficient and the search process will struggle to find the optimum design.
Samareh, J. A., “A survey of Shape Parameterization Techniques”, Tech. rep., NASA Langely Research Center, Hampton, Va., June 1999 (Samareh) presents a survey of possible shape parameterisation techniques, which discusses in detail the advantages and disadvantages of each process. One popular parameterisation technique is to use splines to define the shape. However, for complex geometries this method typically requires a large number of design variables. Free Form Deformation (FFD) is highlighted as a potentially effective shape parameterisation technique. FFD has been successfully implemented as a shape parameterisation method for aerodynamic surfaces within an optimisation process. In particular, FFD is typically used in step 6 of the optimisation strategy described in FIG. 1 in order to alter non-optimal designs. In such a typical optimisation process, FFD is performed directly on the CFD analysis mesh and a CAD geometry is not present within the optimisation loop. If a CAD model of the resulting geometry is required, reverse engineering techniques need to be utilised at step 7 to recreate the CAD model from the deformed CFD mesh. However, it is difficult to recreate the CAD model from the deformed CFD mesh precisely.
Free-Form Deformation (FFD) Techniques and Utilisation
FFD was first introduced in Sederberg, T. W. and Parry, S. R., “Free-Form Deformation of Solid Geometric Models”, SIGGRAPH'86, Vol. 20, ACM, August 1986, pp. 151-160 (Sederberg and Parry) in the field of computer graphics as a method for deforming geometric models in a flexible manner. In general, free-form deformation techniques encompass the following processes: first a volume definition is required and positioned such that it encloses the geometry; a method is then used to deform said volume and finally the geometry inside the volume must deform as part of the volume that encloses it. Sederberg and Parry's FFD method achieves this by using a tricubic Bezier volume that is defined using a lattice of control points. Movements in these control points provide the deformation of the control volume needed. This deformation occurs in an analogous fashion to deforming a Bezier curve through its polygon points. Deformations applied to the Bezier volume are passed to the geometry through its parametric co-ordinates. The following physical analogy, illustrated by FIGS. 2a and 2b, was proposed by Sederberg and Parry and is useful to understand the general FFD method. The geometry to be deformed is in this example a flexible rubber ball 10. This geometry is embedded within a clear flexible plastic block 12. As the block 12 is deformed to produce a deformed block 12a (see FIG. 2b) the embedded geometry 10 is forced to deform in a similar fashion to produce a deformed geometry 10a. The deformed block/geometry 12a, 10a is shown in FIG. 2b. FFD is the mathematical equivalent of this process.
The Sederberg and Parry FFD method can be described in three steps: Control volume/lattice construction; freezing the lattice; lattice and geometry deformation. The three general steps are now discussed in turn.
1. Control Volume/Lattice Construction
As previously mentioned Sederberg and Parry define the control volume which is defined by a regular uniformly subdivided 3D lattice of control points. The lattice consists of l, m and n control points in the S, T and U directions respectively. Within the control volume a local co-ordinate system is defined such that any point X within the control volume is defined byX=X0+sS+tT+uU;  (1)where X0 is the origin of the local co-ordinate system in Cartesian co-ordinates. Note that for any point within the control volume 0<s<1, 0<t<1 and 0<u<1. Although Bezier volumes were used in the discussion by Sederberg and Parry, it is also possible to use other basis functions such as B-splines or non-uniform rational B-splines (NURBS).
Sederberg and Parry restrict the initial lattice shape to be parallelepiped and aligned with the global co-ordinate system. While this limits the types of deformations achievable it has two major advantages. First, because of the convex hull property of Bezier volumes, the parallelepiped condition is necessary to ensure the control volume completely fills the defining lattice. The implication of this is that any geometry enclosed by the lattice can be defined using the Bezier volume parametric co-ordinates. If the initial lattice shape is not parallelepiped the user would need a good knowledge of Bezier curves and their behaviour to form the control volume. The second advantage, which is as a consequence of the first is that the “lattice freezing” process is drastically simplified.
2. Freezing the Lattice
‘Freezing the lattice’ requires the calculation of the position of each geometry point in the local co-ordinate system. For a parallelepiped control volume this can be accomplished through simple linear algebra.
3. Lattice and Geometry Deformation
Deformations are applied to the lattice by moving the control points away from their original positions. The geometry deformed position is then defined by a trivariate tensor product Bernstein polynomial, as
                                          X            ffd                    =                                    ∑                              i                =                0                            l                        ⁢                                          ∑                                  j                  =                  0                                m                            ⁢                                                ∑                                      k                    =                    0                                    n                                ⁢                                                      Bi                    ⁡                                          (                      s                      )                                                        ⁢                                      Bj                    ⁡                                          (                      t                      )                                                        ⁢                                      Bk                    ⁡                                          (                      u                      )                                                        ⁢                  Pi                                                                    ,        j        ,        k                            (        2        )            
In the case of the control volume definition using a Bezier volume the general form of the Bernstein polynomials is:
                                          B            i                    ⁡                      (            s            )                          =                              (                          l              i                        )                    ⁢                                    (                              1                -                s                            )                                      l              -              1                                ⁢                      s            l                                              (        3        )            FFD Advances
Although the FFD method presented by Sederberg and Parry is a powerful deformation tool it does have limitations. The two most significant of these are the restrictions on the lattice shape, and the indirect nature of the control mechanism. Sederberg and Parry restrict the initial lattice shape to be regular, parallelepiped and aligned with the global co-ordinate system. These restrictions limit the types of deformations achievable. The indirect nature of the FFD control mechanism means it is not entirely obvious to the user what deformation can be achieved by moving control points. With respect to shape optimisation this increases the complexity of choosing the design variables needed to achieve the required deformation types. This is further complicated if the control lattice contains internal control points. These limitations have been addressed by numerous studies some of which are discussed below.
1. Basis Functions
The flexibility of the initial control volume representation was increased first in Griessmair, J. and Purgathofer, W., “Deformations of solids with trivariate B-spline”, Eurographics, 1998 (Griessmair and Purgathofer) and then Lamousin, H. J. and Jr, W. N. W., “NURBS-based Free-Form Deformations,” IEEE Computer Graphics and Applications, November 1994, pp. 59-65 (Lamousin and Waggenspack), who extended the FFD method to use B-spline and NURBS bases, respectively. Lamousin and Waggenspack illustrate their method by deforming a model of a human leg. Lamousin and Waggenspack show that the use of NURBS as the basis function enables the use of a non-uniformly spaced lattice and deformation through control point weights. These abilities increase the flexibility of the original FFD method. However, the control lattice is still restricted to a parallelepiped shape.
When using a NURBS basis with the FFD methods the lattice freezing process is complicated as a numerical search technique is needed. This can be avoided if the NURBS weights are all uniform and set to one. Lamousin and Waggenspack use a golden section search routine, whereas a Newton approximation is utilised in Coquillart, S., “Extended Free-Form Deformation: A Sculpturing Tool for 3D Geometric Modeling”, SIGGRAPH'90, Vol. 24, ACM, August 1990, pp. 187-196 (Coquillart). Whichever search process is used can be very computationally expensive. However, if the grid is orthogonal and aligned with the global Cartesian co-ordinates the problem can be separated into three smaller independent searches, one for each of the local co-ordinate directions.
It is important to note that typically during an optimisation process the lattice creation and freezing processes only occur once. New designs are generally created by deforming a ‘base’ design and lattice. Therefore the computational expense required to freeze the lattice is only encountered once.
2. Extended Free Form Deformation (EFFD)
EFFD was first introduced by Coquillart (see above) and expands the FFD methods to use non-parallelepiped shaped lattices. Coquillart uses a tensor product piecewise tricubic Bezier volume to define the geometry deformation. Within each subdivision of a standard l×m×n FFD lattice, termed a “chunk” in Clark, J. H., “Parametric curves, surfaces and volumes in computer graphics and computer-aided geometric design”, Tech. rep., Stanford, Calif., USA, 19 (Clark), Coquillart places an arbitrary tricubic Bezier volume. This extends the original lattice to a (3l+1)×(3m+1)×(3n+1) lattice. Each chunk is equivalent to a single basic FFD volume as presented by Sederberg and Parry. Within Coquillart's EFFD implementation only the corner points of the chunks can be deformed. The positions of the intermediate lattice points are automatically calculated from the chunk point positions and are positioned so as to maintain continuity between the chunk volumes. This is done by maintaining a constant tangent over the chunk boundaries. Freezing the lattice requires first locating the lattice chunk that the geometry resides within and then using numerical search techniques to find the local point co-ordinates. Each geometry point within a chunk is defined by:
                                          X            ffd                    =                                    ∑                              i                ,                j                ,                                  k                  =                  0                                            3                        ⁢                                          Bi                ⁡                                  (                  s                  )                                            ⁢                              Bj                ⁡                                  (                  t                  )                                            ⁢                              Bk                ⁡                                  (                  u                  )                                            ⁢              Pi                                      ,        j        ,        k                            (        4        )            where in this case Bi,Bj,Bk are the control points of fourth order Bernstein polynomials. If the chunk is not parallelepiped a numerical search technique is utilised to solve (4) and find the points.
Coquillart states that the shape of the FFD lattice has a large effect on the deformations achievable. However, the use of lattices which are too complex can produce unpredictable results. The implication of this with respect to shape optimisation is that the choice of the lattice shape is critical to the shape of the design space and consequently the result and efficiency of the optimisation process. A knowledge of the required type of deformation and/or approximate optimum shape will aid lattice selection.
While the EFFD method does enable a greater range of deformations when compared to the original FFD methods, it loses some of the flexibility due to the automatic calculation of intermediate lattice points. It is also clear that the numerical search techniques needed to freeze the EFFD lattice are time consuming and computationally expensive. The EFFD lattice by its very nature will rarely be aligned with the global Cartesian co-ordinate system. Due to this the numerical search cannot be reduced to three line searches, which can be used when a NURBS basis is utilised. Instead a more complicated three variable volume search is required. However, during a design optimisation process freezing the lattice still only needs to occur once to freeze the base design. All subsequent designs will use the frozen base co-ordinates.
3. Directly Manipulated Free Form Deformation (DMFFD)
DMFFD was first introduced in Hsu, W. M., Hughes, J. F., and Kaufman, H., “Direct Manipulation of Free-Form Deformations”, Computer Graphics, Vol. 26, No. 2, July 1992, pp. 177-184 (Hsu et al) as a method of providing a more intuitive interface for the FFD routines. In order to do this a new type of point is defined. An object point is a point placed within the global co-ordinate system on to a particular part of the geometry that can be used as a manipulation handle. The user defines the required position of the object point after deformation. The deformed control point locations needed to produce the required deformation are then calculated and are subsequently used with Eq(3) above to determine the deformed geometry positions. In order to calculate the required control point locations the following systems of equations have to be solved:ΔOb=B*ΔP  (5)
Within the literature there are two distinct methods for solving this equation for the control point displacements. Hsu et al utilised a pseudoinverse of B* enabling the above equation to be rearranged to produce:ΔP=B+ΔOb  (6)
If the number of object points is less than or equal to the number of lattice control points this equation is solvable. However, if the number of object points is greater than the number of lattice control points the system of equations is over-determined and hence there is no unique solution. If the system of equations is under-determined, i.e. the number of control lattice points does not equal the number of object points, multiple solutions exist. In this situation the pseudoinverse gives the least squares solution and so defines the solutions with the minimum control point displacement. In Hu, S.-M., Zhang, H., Tai, C.-L., and Sun, J.-G., “Direct Manipulation of FFD: Efficient Explicit Solutions and Decomposible Multiple Point Constraints”, Visual Computer, Vol. 17, July 2001, pp. 370-379 (Hu et al), a method that uses constrained optimisation is outlined to derive an explicit solution without the use of a pseudoinverse. The resulting solution is:
                              Δ          ⁢                                          ⁢          Pi                ,        j        ,                  k          =                                                    Bi                ,                j                ,                                  k                  ⁡                                      (                                                                  s                        l                                            ,                                              t                        l                                            ,                                              u                        l                                                              )                                                                                                ∑                                      i                    ,                    j                    ,                                          k                      =                      0                                                                            1                    ,                    m                    ,                    n                                                  ⁢                                                      B                                          i                      ,                      j                      ,                      k                                        2                                    ⁡                                      (                                                                  s                        ⁢                                                                                                  ⁢                        1                                            ,                                              t                        ⁢                                                                                                  ⁢                        1                                            ,                                              u                        ⁢                                                                                                  ⁢                        1                                                              )                                                                        ⁢                          (                                                Ob                  t                                -                                  Ob                  1                                            )                                                          (        7        )            
This equation calculates the control point displacements needed to achieve one objective point deformation. Hu et al showed that for multiple objective points the required deformation can be achieved by a series of sequential single objective point calculations, as long as the objective points are mutually independent, i.e., only one object point is present in an individual lattice “chunk”.
With respect to shape optimisation, when DMFFD is used the design variables are the objective point deformed co-ordinates. This method can be used to reduce the number of design variables while maintaining the flexibility within the geometry representation.